scripts/Empirical example NatMath/NatMathOptimalBin.R

In joakimwallmark/PolyOptimalIRT:

### Only Binary items!

# ----------- preliminaries -----------
library(TestGardener)

# ----------- read in data -----------
titlestr  <- "National Math Test"
U         <- scan("data/NatMath.txt","o")
N         <- length(U) # Number of examinees
Umat      <- as.integer(unlist(stringr::str_split(U,"")))
n         <- length(Umat)/N # Number of items
U         <- matrix(Umat,N,n,byrow=TRUE)

#  data in score mode, convert to index mode
U = U + 1

key     <- rep(1, 12)

# Define the option score values for each item ----------------------------

binindices <- c(1:5, 7:9, 11:13, 16)
Ubin <- U[, binindices]

optList <- vector("list",32) # option scores
optList[[ 1]] = c(0, 1)
optList[[ 2]] = c(0, 1)
optList[[ 3]] = c(0, 1)
optList[[ 4]] = c(0, 1)
optList[[ 5]] = c(0, 1)
optList[[ 6]] = c(0, 1)
optList[[ 7]] = c(0, 1)
optList[[ 8]] = c(0, 1)
optList[[ 9]] = c(0, 1)
optList[[10]] = c(0, 1)
optList[[11]] = c(0, 1)
optList[[12]] = c(0, 1)

itemLab <- vector("list", 32)
itemLab[[ 1]] <-  'Question  1 '
itemLab[[ 2]] <-  'Question  2a'
itemLab[[ 3]] <-  'Question  2b'
itemLab[[ 4]] <-  'Question  3a'
itemLab[[ 5]] <-  'Question  3b'
itemLab[[ 6]] <-  'Question  5 '
itemLab[[ 7]] <-  'Question  6 '
itemLab[[ 8]] <-  'Question  7 '
itemLab[[9]] <-  'Question  8b'
itemLab[[10]] <-  'Question  9a'
itemLab[[11]] <-  'Question  9b'
itemLab[[12]] <-  'Question 12a'

optListBin <- list(itemLab=itemLab, optLab=NULL, optScr=optList)

NatMath_dataList <- make.dataList(Ubin, key, optListBin, scrrng = c(0,12))

#  Set initial values that are required in the later analysis -----------------
#  compute the initial option surprisal curves using the
#  percentage ranks as initial estimates of theta
theta     <- NatMath_dataList$percntrnk
thetaQnt  <- NatMath_dataList$thetaQnt
# thetaQnt  <- NatMath_dataList$thetaQnt[rep(c(T, F), 23)]
chartList <- NatMath_dataList$chartList
WfdResult <- Wbinsmth(theta, NatMath_dataList, thetaQnt, chartList)

#  Plot the initial option proability and surprisal curves

WfdListinit <- WfdResult$WfdList
binctrinit  <- WfdResult$aves
Qvecinit    <- c(5,25,50,75,95)

Wbinsmth.plot(binctrinit, Qvecinit, WfdListinit, NatMath_dataList, Wrng=c(0,7))
Wbinsmth.plot(binctrinit, Qvecinit, WfdListinit, NatMath_dataList, Wrng=c(0,7),
              plotindex=20)

#  Proceed through the cycles---------------------------------------------------
ncycle=10
AnalyzeResult <- Analyze(theta, thetaQnt, NatMath_dataList, ncycle=ncycle, itdisp=FALSE)

#save(NatMath_dataList, AnalyzeResult, file = "data/NatMath_fittedmodel_bin.RData")
parList  <- AnalyzeResult$parList
meanHvec <- AnalyzeResult$meanHvec

#  ----------------------------------------------------------------------------
#              Plot meanHsave and choose cycle for plotting
#  ----------------------------------------------------------------------------

cycleno <- 1:ncycle
par(mfrow=c(1,1))
plot(cycleno,meanHvec, type="b", lwd=2, xlab="Cycle Number")

#  select cycle for plotting

icycle <- 10

NatMath_parListi  <- parList[[icycle]]

WfdList    <- NatMath_parListi$WfdList
Qvec       <- NatMath_parListi$Qvec
binctr     <- NatMath_parListi$binctr
theta      <- NatMath_parListi$theta
arclength  <- NatMath_parListi$arclength
alfine     <- NatMath_parListi$alfine

#  ----------------------------------------------------------------------------
#                   Plot surprisal curves for each test question
#  ----------------------------------------------------------------------------

#  plot both the probability and surprisal curves along with data points
Wbinsmth.plot(binctr, Qvec, WfdList, NatMath_dataList, Wrng=c(0,7), saveplot = F)

#  ----------------------------------------------------------------------------
#                         Plot density of theta
#  ----------------------------------------------------------------------------

ttllab     <- paste(titlestr,": percent rank", sep="")
edges      <- c(0,100)
theta_in   <- theta[0 < theta & theta < 100]
indden10   <- scoreDensity(theta_in, edges, 15, ttlstr=ttllab)